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ABSTRACT 

Using wide field narrow-band surveys, we provide a new measurement of the z = 6.6 Lyman- 
a Emitter (LAE) luminosity function (LE), which constraints the bright end for the first time. 
We use a combination of archival narrow-band NB921 data in UDS and new NB921 mea¬ 
surements in SA22 and COSMOS/UltraVISTA, all observed with the Subaru telescope, with 
a total area of ^ 5 deg^. We exclude lower redshift interlopers by using broad-band optical 
and near-infrared photometry and also exclude three supernovae with data split over multiple 
epochs. Combining the UDS and COSMOS samples we find no evolution of the bright end 
of the Lya LE between z = 5.7 and 6.6, which is supported by spectroscopic follow-up, 
and conclude that sources with Himiko-like luminosity are not as rare as previously thought, 
with number densities of ~ 1.5 x 10“® Mpc“^. Combined with our wide-held SA22 measure¬ 
ments, our results indicate a non-Schechter-like bright end of the LE at z = 6.6 and a different 
evolution of observed faint and bright LAEs, overcoming cosmic variance. This differential 
evolution is also seen in the spectroscopic follow-up of UV selected galaxies and is now also 
conhrmed for Lya emitters, and we argue that it may be an effect of re-ionisation. Using a toy- 
model, we show that such differential evolution of the LE is expected, since brighter sources 
are able to ionise their surroundings earlier, such that Lya photons are able to escape. Our tar¬ 
gets are excellent candidates for detailed follow-up studies and provide the possibility to give 
a unique view on the earliest stages in the formation of galaxies and re-ionisation process. 

Key words: galaxies; high-redshift, galaxies: luminosity function, cosmology; observations, 
galaxies: evolution, cosmology: dark ages, re-ionisation, hrst stars. 


1 INTRODUCTION 

The Lyman-a (Lya) emission line (1216 A) is a powerful tool to 
study the formation of galaxies in the early Universe. This is be¬ 
cause it has been predicted to be emitted by young “primeval” 
galaxies jPartridge & Peebles||1967{ |Pritchet||l994| ), but also be¬ 
cause it is redshifted into optical wavelengths at z > 2, where most 
rest-frame optical emission lines are impossible to observe with 
current instrumentation. 


Indeed, the Lya line has been used to spectroscopically con¬ 
firm high redshift candidate galaxies up to z ~ 7.5 obtained 


with the Lyman-break technique (e.g. 

Steidel et al.|1996| |Finkel-| 

Istein et al.|2013||Schenker et al.|2014 

1 , which is based on broad- 


band photometry using e.g. WFC3 on the Hubble Space Telescope 
(HST). Galaxies selected this way are called Lyman-Break Galax- 
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ies (LBGs) and the current largest sample contains already 10,000s 
(e.g. |Bouwens et al.|2014b| >. 

Narrow-band surveys select emission line-galaxies at specific 
redshift slices and are therefore used to search for Lya emitters 
(LAEs) directly. Samples of LAEs have now been established from 
z ~ 2 — 7 through narrow-band surveys (e. 


Rhoads et al.|2000 

IFynbo et al.|2001||Rhoads et al.|2003| 

Vlalho- 

tra & Rhoads|2004 

[|Taniguchi et al.|2005| Shimasaku et a 

2006 

Westra et al.|2006 

INilsson et al.|2007l lOuchi et al.||2008 

2010 

Hu et al.|2010||Kashikawa et al.|2011[|Shibuya et al.|2012| 

Konno 


VUDS and MUSE; [Hiiret al.|2008||Cassata et al.|20l5l|Bacon et al.| 
|2014[ l. Limited samples of LAEs at lower redshifts and the local 
Universe also exist that are detected through e.g. GALEX or HST 
(e.g. |Hayes et al.|2007[|Deharveng et al.|2008[[Cowie et al.|2010^ . 

While part of the difference between LAEs and LBGs is just 
the way they are detected, there are also differences in their av¬ 
erage properties. There exists an anti-correlation between the UV 
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brightness and the Lya Equivalent Width (EW) ( |Ando et al.|2006 [ 
|Stark et al.||2010) , indicating that the brightest LBGs are typi¬ 
cally not LAEs, and that the UV continuum for most LAEs is 
very hard to detect, even in the deepest broad-band images ( |Bacon| 
|et al.|2014) . Spectroscopic follow-up of LBG selected galaxies has 
shown that the typical Lya EW increases with increasing redshift 
up to 2 ~ 6.5. This is likely due to LBGs being younger on aver¬ 
age and less dustier at higher redshift (jStark et al.|2010| |Schenker| 
|et al.|20i2{|Cassata et al.|2Q15^ . This picture is consistent with the 
evolution of the luminosity function (LF) of the different classes 
of galaxies. For LAEs, the Lya LF is remarkably constant between 
a = 3 — 6 (e.g.|Shimasaku et al.|2006||Dawson et al.|2007||Gron-| 
[wall et al.||2007[ |Ouchi et al.||2008t , while the UV LF of LBGs 


declines to higher redshifts in a reasonably uniform way due to the 
decline in the global star formation activity in galaxies jEllis et ^ 
|2013||Bouwens et al.|2014b}[McLeod et al.|2014t . This also indi¬ 
cates that the Lya emission line generally brightens with increasing 
redshift. 

From these observations, the picture has emerged that Lya is 
preferentially observed at a specific phase during a galaxy’s evolu¬ 
tion. Since Lya is produced by recombination radiation from hy¬ 
drogen clouds around very massive, young (< 10 Myr) stars (e.g. 
|Schaerer|2003^ , and Lya is easily absorbed and re-scattered (lead¬ 
ing to lower surface brightnesses), on average LAEs are believed 
to be young starbursts, while LBGs on average are slightly more 
evolved galaxies with a higher dust content (e.g. |Verhamme et al.| 
|2008|l.|0no et al.|(|2010|l find that the UV slope of z = 6 — 7 LAEs 
is very steep (/( = —3), while [Bouwens et al.| ( |2014a) find that 
the UV slope of LBGs at similar redshifts is typically slightly shal¬ 
lower (/3 = —2.2). From clustering measurements, [Gawiser et aT] 
( |2007^ ; [Ouchi et ar]j2010| l and |Bielby et al.H2015 l agree on an 
average LAE halo mass of ~ 10^^ M© from 2 ; = 3 — 7. For LBGs 
alternatively, the typical halo mass is one order of magnitude higher 
at these redshifts (e.g.|Ouchi et al.|2003[[Ha 


ietal.|2004[[0uchi| 


let al.|2005HHildebrandt et al.|2009t , more typical of “Milky Way” 
dark matter haloes of 10^^ Mq. 

Near the re-ionisation redshift, physical processes start to play 
a role which are additional to intrinsic changes in the properties of 
galaxies, since Lya is easily absorbed by a neutral Inter Galactic 
Medium (IGM). While LAEs can be an important source of ionis¬ 
ing photons for re-ionisation, one of the main interests in studying 
LAEs at this epoch is observable effects of a higher neutral IGM 
opacity. Besides evolution of the luminosity function, these ob¬ 
servables include an increased observed clustering in a more neu¬ 
tral IGM, and attenuated line-profile. The observed clustering of 
LAEs increases since the observability is favoured when sources 
are in overlapping ionised spheres ( |McQuinn et al.|[2007| ), while 
the line-profile becomes more asymmetric due to absorption and 
re-scattering of Lya photons in a more neutral medium dDijkstral 
let al.|2007| l. 


At 2 > 7 spectroscopic follow-up of LBGs is remarkably less 


successful IFontana et al.|2010l 

Stark et al.|2010 Pentericci et al. 

201 mono et al.|2012l|Treu et al. 

2013l|Pentericci et al.|2014l|Caru- 

ana et al.|2014J, indicating eithc 

;r a lower intrinsic escape of Lya 
reased column density of absorb- 
2013 1 , or a higher neutral fraction 
a et al.|2007||Schenker et al.|2014| 
2014J. There is also evidence for 

(e.g. Dijkstra et al.|2014J, an inc 
ing clouds (Bolton & Haehnelt| 
of the IGM I Santos|2004| Dijksti 
Taylor & Lidz|2014| Tilvi et al. 


an increased opacity to Lya photons from the Lya LF, which is 


observed to decline very rapidly ( |Ouchi et al.|2010{ |Konno et al.| 
|2014| >. Searches for LAEs at 2 = 7.7 and 2 = 8.8 have been 
unsuccessful in spectroscopically confirming any of the candidates 


(Willis & Courbin|20051 Cuby et al.|2007| Willis et al.|2008||So- 

bral et al.|2009l|Hibon et al.|2010 

ITilvi et al.|2010||Clement et al. 

2012 

|Krug et al.|2012[|Jiang et al. 

20131 Faisst et al.|2014| Matthee 

et al. 

2014|(. However, these studies are still limited by their sensi- 


tivity since Lya is shifted into the near-infrared (NIR). Most of 
these studies only probe tiny areas in the sky, meaning that bright 
sources might be missed. 

Typically, research is so far limited to ~ 1 deg^ areas (e.g. 
|Ouchi et al.||200^ |2010| >, where cosmic variance, especially for 
the observability of Lya around the re-ionisation epoch, can play 
a large role. To make further progress, we are carrying out an 
extensive set of wide-field narrow-band observations to study 
the evolution of Lya emitters from the epoch of re-ionisation 
(2 ~ 6 — 9) up to the peak of the cosmic star formation history 
(2 ~ 2 — 3). Our aim is to explore the evolution of the bright end 
which is so far uncharted and for which spectroscopic follow-up 
is easier and gives a better comparison to to surveys at the highest 
redshifts. 


In this paper, we focus at the 2 = 6.6 Lya LF because of its 
importance to the study of re-ionisation. The widest narrow-band 
survey at that redshift to date has been presented by |Ouchi et al.| 
(2010 1 , which reaches a Lya luminosity of ~ 10’’^'®erg over a 
~ 0.9 deg^ area. The brightest Lya emitt er in their sample, Himiko , 
with a luminosity of 3.5 x 10'*^ erg s~^l Ouchi et al.|2009 2013 1 , 
has been seen as a very rare source, a triple merger, one of its kind. 
Because there has been only one very bright source known, the er¬ 
ror on its number density is large, such that there is a factor 30 
offset between the fitted luminosity function and the data. We have 
obtained wide field observations to further constrain the number 
density of bright Lya emitters. We introduce our data and sam¬ 
ple selection in §2. Using our data set, we re-produce the |Ouchi| 
|et al.| ( [2()10^ sample and add new 2 = 6.6 Lya candidates using 
deep archival Subaru data in the COSMOS field and new shallow 
wide-field data in SA22 in §3. Our estimates of the completeness 
of our selection procedure and description of the corrections made 
to the luminosity function are shown in §4. This leads to a new es¬ 
timate of the luminosity function in §5 (supported by spectroscopic 
confirmation of the two brightest LAEs in COSMOS ; |Sobral et al.| 
( |2015b| l), where we find that the combined luminosity function has 
a non Schechter-like bright end. We discuss the evolution and im¬ 
plication for re-ionisation in §6. The main results are shown in Fig. 

S which shows our new estimate of the 2 = 6.6 Lya LF, and Fig. 

where we compare the evolution between 2 = 5.7 and 2 = 6.6. 
Throughout the paper, we use a 737 ACDM cosmology (Hq = 70 
km s“^Mpc“^, Om = 0.3, Oa = 0.7) and magnitudes are mea¬ 
sured in 2” diameter apertures in the AB system. 


2 OBSERVATIONS & DATA REDUCTION 
2.1 Imaging 

Optical imaging data was obtained with Subaru’s Suprime-Cam 
( |Miyazaki et al.|2002) using the NB921 narrow-band filter (Fig.[^. 
The Suprime-Cam is composed by 10 CCDs with a combined field 
of view of 34' x 27' and with chip gaps of ~ 15". The NB921 
filter has a central wavelength of 9196 A and a FWHM of 132 A 
and is located in a wavelength region free of OH line-emission in 
the atmosphere (see Fig.[^. 

We obtained archival ultra-deep observations in UKIDSS- 
UDS (02 18 00 -05 00 00) and COSMOS-UltraVISTA (10 00 00 











































































































The Lya luminosity function a? z = 6.6 3 


Field 

R.A. 

(J2000) 

Dec. 

(J2000) 

Int.time 

(ks) 

FWHM 

(”) 

Area 

(deg2) 

Depth 

(3a) 

Dates 

COSMOS-1 

10 01 28 

- 1-02 25 51 

10.8 

0.3 

0.24 

25.8 

2009 Dec 19 

COSMOS-2 

09 59 35 

-1-02 27 01 

8.78 

0.3 

0.24 

25.9 

2009 Dec 19,20 

COSMOS-3 

10 01 24 

-1-01 58 00 

10.8 

0.6 

0.24 

25.9 

2009 Dec 21 

COSMOS-4 

09 59 29 

-1-01 58 42 

7.80 

0.6 

0.24 

25.7 

2009 Dec 21 

SA22-DEEP 

22 19 14 

-1-00 11 24 

32.1 

0.3 

0.17 

26.4 

2009 Sep 15-17 

SA22-WIDE-11-19]* 

22 16 19 

-1-00 10 00 

0.36 

0.5 

2.72 

24.3 

2014 May 28-29 

UKIDSS-UDS C 

02 18 00 

-05 00 00 

30.0 

0.8 

0.14 

26.4 

2005 Oct 29, Nov 1, 2007 Oct 11,12 

UKIDSS-UDS N 

02 18 00 

-04 35 00 

37.8 

0.9 

0.18 

26.4 

2005 Oct 30,31, Nov 1, 2006 Nov 18, 2007 Oct 11,12 

UKIDSS-UDS S 

02 18 00 

-05 25 00 

37.1 

0.8 

0.19 

26.4 

2005 Aug 29, Oct 29, 2006 Nov 18, 2007 Oct 12 

UKIDSS-UDS E 

02 19 47 

-05 00 00 

29.3 

0.8 

0.16 

26.4 

2005 Oct 31, Nov 1, 2006 Nov 18, 2007 Oct 11,12 

UKIDSS-UDS W 

02 16 13 

-05 00 00 

28.1 

0.8 

0.14 

26.4 

2006 Nov 18, 2007 Oct 11,12 


Table 1. Observation log for the NB921 optical imaging in the COSMOS, SA22 and UKIDSS UDS fields. Depths are based on empty aperture measurements 
and take correlations in the background into account. The UDS data have been analysed by |Ouchi et al.|{^10^ . The area is the area after masking (see §2.1). 
* Note that we have 19 pointings in SA22-WIDE, all with similar observing conditions. 
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Figure 1. Filter transmission profile of the NB921(blue) and z (green). The 
transmission is normalised to the maximum transmission in each filter. In 
grey, OH emission lines from the night sky are shown pousselot et al.| 
|2000) . It can be seen that the naiTow-band is in sky-line free region allowing 
us to go very deep. This also facilitates the spectroscopic follow-up. 

-1-02 00 00) and SA22 (22 00 00 -l-OO 00 00) and took new data in 
a wide area (~ 4.5 deg^) in SA22 on May 28-29 2014, observing 
program S14A-086 (PI: Sobral). These fields were chosen for their 
multi-wavelength coverage and low Galactic foreground emission. 
Apart from providing a large total area, the combination of ultra¬ 
deep and shallower surveys allows us to sample a wide range of 
luminosities. 

In UDS, there are 5 sub-fields with a total integration time ranging 
from 7.8 to 10.5 hours (see Table 1). This exposure time is obtained 
after stacking individual exposures of 1.2ks, with a small dithering 
pattern. The seeing FWHM ranges from 0.8 — 0.9”. We mask re¬ 
gions around bright stars, where spherical artefacts boost the fluxes 
artificially. We also mask horizontal and vertical stripes caused by 
blooming of a saturated bright object (see Fig.[^. This is the same 
raw data as used to study Lya emitters at 2 : = 6.6 by|Ouchi et al.| 

d^Ot . 

In the four pointings in COSMOS, the total integration time ranges 
from 2.2 to 3 hours, with individual exposure time of 1.2ks, such 
that the total number of exposures per field is smaller (ranging from 
7 to 9). In two pointings, the seeing is particularly good (0.3”), the 
other two have seeing FWHM of 0.6”. Similar to in UDS, we mask 
spherical halo-regions around bright stars, and vertical and hori¬ 
zontal blooming-stripes (depending on the position angle rotation 
of the Suprime-cam pointing). 


The deep SA22 data consists of 27 exposures of 1.2ks at a single 
pointing with small dithering in perfect seeing conditions (FWHM 
0.3”). We mask a small noisy region and blooming patterns. The 
wide SA22 data consists of 19 pointings of 3 exposures of 120s 
in very good seeing conditions (FWHM 0.5”). Because of the lim¬ 
ited number of exposures, a significant area has a lower signal to 
noise due to the dithering pattern. Another minor issue is that the 
astrometric corrections and calibration of the zero-point are less 
accurate in certain detectors due to the low signal to noise. We con¬ 
servatively mask all these regions and also mask regions around 
bright stars and are left with a final coverage of 2.7 deg^. 

The total area after masking is 4.66 deg^. An overview of the ob¬ 
servations is given in Table 1. 

2.2 Data reduction 

The NB921 imaging data was reduced with SDFRED2 package 
( |Ouchi et al.|2004| . The sequential steps in the reduction are: 

(i) Overscan and bias subtraction: for each image, a median 
value for the overscan region was determined and subtracted in 
each line of pixels. The bias was subtracted by assuming that it 
has the same value as the overscan. 

(ii) Flat fielding: flat frames are obtained by observing a uni¬ 
form light source and are required to correct variations in pixel- 
to-pixel sensitivity across the camera. By dividing the images by 
these flat frames, the background becomes flat and luminous pat¬ 
terns caused by differences in the sensibility of the pixels are re¬ 
moved. 

(iii) Point spread function homogenisation: the point spread 
function (PSF) measures the response of the detector to a point¬ 
like source. PSF sizes were obtained by measuring the FWHM of 
the point-like sources of each frame. The target PSF for homogeni¬ 
sation was defined as the one with more occurrences and the frames 
with PSF smaller than the target were smoothed with a gaussian. 

(iv) Sky background subtraction: a mesh pattern was com¬ 
puted to represent the sky background, the pattern was interpolated 
and subtracted on each frame. 

(v) Bad pixel masking: defects with the detector and problems 
with the observations may cause data in some pixels to become 
corrupted. A mask is applied to these pixels. 

(vi) Astrometric calibrations: we correct each image for astro- 
metric distortions using Scamp ( |Bertin|2006| l, which fits a polyno¬ 
mial solution by matching detected sources with the 2MASS cata¬ 
log in the J band ([Skrutskie et al.|2006^. It also takes into account 
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that images have different integration times and attributes a differ¬ 
ent weight to each image. 

(vii) Stacking: Once each individual frame has been reduced, 
we stack the different jittered frames for each pointing. 

(viii) Cosmic ray rejection: cosmic rays are rejected automat¬ 
ically based on the standard deviation in the pixel values in a 1” 
aperture. The standard deviation is typically a hundred times higher 
for cosmic rays than for real sources. We use a very conservative cut 
since we do not want to risk rejecting real sources, meaning that our 
sample will still be somewhat contaminated by cosmic rays. Partly 
because of this, we inspect all our final candidates visually. 


2.3 Photometric calibration & survey depths 

Once we obtained the reduced data for each pointing, we set 
the zero-point (ZP) to a magnitude of 30 (AB). This is done by 
extracting sources with SEXTRACTOR jBertin & Amouts||1996^ 
with high detection thresholds (> 20a) and match these sources 
with public catalogues in UDS ( [Cirasuolo et al.|2007| l, COSMOS 
pibert et al.||2009T l and our own catalog in SA22 (based on K 
detected sources in UKIDSS-DXS; [Matthee et al.|[^14[ jSobralj 
jet al.|[^15aj , using STILTS ( [Taylorj [2006^ . We only include 
sources with narrow-band magnitudes brighter than 19, such that 
our detections are at sufficient high signal to noise, and fainter 
than 16, since brighter sources are saturated in our data. In each 
pointing, we use roughly 500 sources. We then set the ZP by 
correcting the cropped mean difference between the magnitudes in 
our images and the ones in the catalogue. 

We estimate our survey depth by measuring the root mean 
square (rms) of the background in a million empty apertures with 
2” diameters, placed at random places in the image, avoiding 
sources which are detected at > 3a. Empty aperture measurements 
take into account that the background noise is correlated and are 
thus more robust than if the background is measured on a pixel by 
pixel basis (c.f. jMilvang- Jensen et al.|2013| l. Also, empty apertures 
can still include very faint sources below our detection threshold, 
so this is a conservative upper limit. The depths of the narrow-band 
images are listed in Table 1. The 3a depths are 26.4 in UDS, 25.8 in 
COSMOS, 26.4 in the deep pointing in SA22 and 24.3 in the wide 
pointings in SA22. 

2.4 Multi-wavelength data and photometry 

Eor UDS, deep z-band data (26.6, 3a) is available from the Subaru 
Extreme Deep Eield (SXDE) project ( [Furusawa et al.|[2008[ l, as 
well as data in the optical bands B, V, R and i, with 3a limits: 
28.3, 28.4, 27.8 and 27.2 respectively (see Table 2.). This multi¬ 
wavelength data is essential to identify different line-emitters. 
The images in the optical and NB921 are all aligned, because all 
are from a single survey, telescope and instrument. Furthermore, 
UKIDSS NIR J, H, K data ( [Lawrence et al.|2007[ l is available for 
60 % of the coverage, with 3cr depths 25.4, 24.7 and 24.9. For NIR 
photometry, we use SWARP jBertin|20I0) to align the NIR images 
to the NB92I images and interpolate the NIR images, since the 
pixel scale is slightly larger (UKIRT WFCAM; |Casali et al.|2007^ 
than the Suprime-Cam pixel scale. 


imaging, with 3a depths of 27.6, 27.0, 26.9 and 25.8 (Table 2) 
which is available through the COSMOS archiv^ We align the 
optical images to the narrow-band images using SWARP. NIR data 
in Y JHKs is available from Ultra VISTA DR2 jMcCracken et al.[ 
[20I2| | with 5cr depths 25.4, 25.1, 24.7 and 24.8. The pixel scale 
of VISTA’s VIRCAM is 0.15”/px, so we degrade the images to 
the pixel scale of the narrow-band images (0.2’Vpx) and align them. 

The SA22 field is covered by CFHTLS and UKIDSS DXS 
surveys. SA22 is W4 in CFHTL^and is imaged in ugriz with 
MegaCam, which has a field of view of 1 x 1 deg^ and a pixel 
scale of 0.187”/pixel. The near-infrared UKIDSS DXS survey 
has imaged it in J and K filters with UKIRTAVFCAM. Note 
that the multi-wavelength data is not as deep as in the other two 
fields (typically 1-2 magnitudes shallower), and there is also no 
Spitzer/IRAC data available. This limits the search for LAEs in the 
deep pointing since the uncertainty in the z-band is much higher 
than the uncertainty in NB92L For the brightest objects (including 
all reliable detections in the Wide coverage), this is less of a 
problem. As before, we align the optical and near-infrared images 
to our narrow-band pointings using SWARP, including degrading 
the pixel scale to that of the narrow-band imaging. 

For all fields, photometry is extracted using SEXTRACTOR in 
dual-mode with the NB921 image as detection image and within 
a 2” circular aperture. In the case of a non-detection in any of 
the broad-bands, we assign la limits. We then use the i-band to 
correct our 2 -band such that the median narrow-band excess is 
zero for all sources by fitting a linear relation between the {i — z) 
colour and the narrow-band excess. Since the narrow-band filter 
is almost in the center of the broad-band filter, the correction is 
small, Zcor = z — 0.13(i — z) + 0.286. For sources undetected in 
the i band we assign the median correction of -1-0.03. 


3 SELECTING LYMAN-a EMITTERS AT Z = 6.6 

In the NB92I data, Lyman-a emitters need to be selected as 
line-emitters at 2 = 6.55 ± 0.055. This requires multi-wavelength 
coverage of the fields, which is available through a combination of 
large (public) surveys. First of all, we require broad-band photom¬ 
etry over the same wavelength coverage as the narrow-band, which 
is in this case the 2 -band. 

Line-emitters are selected based on two criteria (e.g. |So-[ 
[bral et ^[2013t : the first is that the narrow-band excess must 
be high enough. Since the observed equivalent width of emission 
lines increases with redshift and the intrinsic EW of Lya is high 
(~ 100 - 200 A), we expect Lyman-a emitters to have a high 
narrow-band excess. We follow previous searches (e.g. [Ouchi et al.[ 
[2010| ( and use an excess criterion of 2 — NB021 > 1, correspond¬ 
ing to a 2 = 6.6 rest-frame EW of 38 A. This limit is also chosen to 
minimise contamination by lower redshift interlopers, although we 
will also lower this criterion and comment on the differences. To 
convert the narrow-band excess to the observed EW, we first trans¬ 
form magnitudes (mU to flux densities in each filter (/U with the 
standard AB convention: 


The COSMOS field is one of the best studied extra-galactic 
fields with > 30 bands coverage jllbert et al.[[200^ , ranging 
from X-ray to radio. We use deep optical BViz data from Subaru 


Capak et al. 1 2007 i; http://irsa.ipac.caltech.edu/data/COSMOS/ 
http://www.ctht.hawaii.edu/Science/CFHTLS/ 
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Field 

Surveys 

Optical 

Depths 

Near — infrared 

Depths 

UDS 

SXDF, UDS & SpUDS 

BVRiz 

28.3, 28.4, 27.8, 27.2, 26.6 

JHK 

25.7, 24.5, 24.4 

COSMOS 

COSMOS, UltraVISTA & S-COSMOS 

BViz 

27.6, 27.0, 26.9, 25.8 

YJHK 

25.9, 25.4, 25.1,24.7 

SA22 

CFHTLS & DXS 

ugriz 

26.2, 26.3, 26.0, 25.7, 24.5 

JK 

24.4, 23.9 


Table 2. 3cr depths of broad-band coverage of our survey-fields obtained using empty aperture measurements. Note that there can be more wavelength data 
available, but we only use the broad-bands listed here for consistency between the helds. We obtain our own photometry by first registering all the images to 
the Subaru NB921 measurements and then extracting photometry within 2” apertures in dual mode. The SXDF data is presented in|Furusawa et al.|j2008|, the 
COSMOS data in |Ilbert et al.H20Q^ and Ultra VISTA in |McCracken et ar]j2012) . 


fi = 


C 1 ^-0.4(mi+48.6) 

\2 

center 


( 1 ) 


In this equation, c is the speed of light and Xi^center is the 
central wavelength in each filter, which are 9183.8 A and 8781.7 A 
for the narrow-band (NB) and broad-band (BB) respectively. 
Using Eq. 1, we use the following equations to convert to EW and 
line-flux respectively: 


EW = AXnb 


/nb —Jbb 

r p A A 

JBB - JNB-^ 


NB 

BB 
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Here, Jmb and fsB are the flux-densities, AXmb and A Ass 
the filter-widths, 135.1 A and 1124.6 A, respectively. In this for¬ 
mula, the numerator is the difference in narrow-band and broad¬ 
band flux and the denominator the continuum, which is corrected 
for the contribution from the narrow-band flux. The formula breaks 
down at a certain NB excess depending on the specific filters, and 
thus we set the EW of those sources to > 1500 A. 

The line-flux is computed using: 


fline = AXnB 


/nB — fBB 

1 _ AAjvg 

-L A \ 


(3) 


The second criterion for selecting emission line galaxies is that 
the excess should be significant, meaning that it is not dominated 
by errors in the narrow-band and broad-band photometry. We will 
follow the methodology presented in |Bunker et al.| (T995| > and the 
equation from |Sobral et al.| ( |2013[ l to compute the excess signifi¬ 
cance (S): 


(NIR) J, H and K filters, can provide valuable insight in the na¬ 
ture of the candidates and possibly help excluding lower-redshift 
interlopers. Finally, Spitzer-IRAC data can be used as further con¬ 
straints on excluding dusty low redshift interlopers (generally with 
bright IRAC detections and red colours), or as a further evidence 
for the source being at 2 = 6.6, since at that redshift the [3.6] and 

[4.5] pm bands are contaminated in such a way that sources with 
strong nebular emission (EW > 1000 A) are expected to have blue 

[3.6] - [4.5] colours (e.g. jStark et al.|2013|[SmIt et al.|2014a|b) . 

In addition to using colours for our selection and characterisation 
of Lya candidates, we also compute photometric redshifts using 
EAZY vl.l jBrammer et al.|[2008t , which includes the contribu¬ 
tion of emission lines (although typically not strong enough for the 
observed extreme emission lines galaxies). We include optical and 
near-infrared photometry and we use it to identify possible lower 
redshift interlopers. 

The details of our selection per field are presented in §3.1 for UDS, 
§3.2 for COSMOS and §3.3 for SA22. 


3.1 Selecting LAEs in UDS 

We use our photometry described in §2.4 and select line-emitters 
using the following criteria: 

(i) z - NB921 > 1 

(ii) S > 3 

(iii) Pass visual inspection 


1 _ 1 r)-OA(BB-NB) 

S = -- , - (4) 

In this case, BB is the z-band magnitude after correction using the 
i band (see next subsections), NB the magnitude in NB921, ZP 
the zero-point of the images, which is set to a 30 AB magnitude. 
apx is the root mean squared (rms) of background pixel values in 
the data of the respective filters and Tap is the aperture radius in 
pixels. Our depths are estimated using empty aperture based rms 
values, which takes correlations in the background into account. 
For the selection of emitters however, we use pixel based rms 
values for consistency with previous surveys, but also check that 
our results are robust when using empty aperture-based E values. 

After selecting a sample of line-emitters, we use multi¬ 
wavelength data to distinguish high redshift candidate LAEs. In 
addition to the z-band, we also need bands in bluer wavelengths 
in order to apply the Lyman break technique to select high redshift 
sources. In this case this means that there should be no detection in 
the B, V, u, g and r filters and a strong break in the (i — z) colours. 
Measurements in redder wavelengths, such as in the near-infrared 


After the first two criteria, we have 1514 line-emitter candidates 
(see Table 3 for the numbers after each step). We check each line- 
emitter candidate visually in the NB921 image and exclude 122. 
Furthermore, 25 candidates are excluded since their excess is un- 
physically high, meaning that their (non-)detection in the 2 -band is 
in disagreement with the flux solely contributed by the measured 
narrow-band flux by more than 3a. Most of the excluded candi¬ 
dates have their flux boosted by artifacts from bright stars, such as 
haloes or spikes (even after masking), others are excluded because 
they reflect read-out noise (since the images are so deep) and there 
are CCD-grid like patterns. 

In total, we find 1367 line-emitters, which selection is shown 
in Fig. 1^ This sample is dominated by H/3/[OIII] at 2 = 0.83, 
[Oil] emitters at 2 = 1.46 and Hq at 2 = 0.40 (based on photo¬ 
metric redshifts), even though the high excess criterion is already 
used to minimise this number. These lower redshift emitters are 
described for example in |Sobral et al.| ( |T013[ l, which also shows a 
photometric redshift distribution. To select Lya candidates, we ex¬ 
clude low redshift sources and select high-redshift candidates using 
the Lyman-break technique: 
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Figure 2. Narrow-band excess versus narrow-band magnitude. This figure shows the selection of line-emitters in our different fields. The grey points show all 
detected objects (after removing cosmic rays, visually identified spurious sources and sources for which the excess is unphysical) and the green points show the 
selected LAEs in each field. The selection consists of an EW cut (solid horizontal line) and a criterion to determine whether the excess is significant (dashed 
lines). The EW cut comesponds to a restframe EW of 38 A at z = 6.56, and the significance cut corresponds to 3S. After visual removal of spurious sources 
and cosmic rays and optical and near-infrared broadband criteria, we find 99 LAEs in UDS, 15 in COSMOS, 2 in SA22-Deep and 19 in SA22-Wide. While 
the narrow-band imaging in UDS and SA22-Deep have similar depth, the number in UDS is higher due to a roughly 5 times larger area and (most importantly) 
deeper z-band imaging. There are less candidates in COSMOS due to shallower narrow-band imaging. The imaging in SA22-Wide is even shallower, but this 
is compensated by a much larger area. 


B > 28.7 A y > 28.2 A (i - z > 1.3 V i > 27.2) 

The Lyman-break technique is based on the absence of flux blue- 
ward of the Lyman limit (912 A) and we therefore require our 
sources not to be detected at wavelengths below 912 x (6.6 + 1) = 
6930 A. This results in non-detections in B and V, and a strong 
{i — z) break (or a non-detection in i as well). The {i — z) cri¬ 
terion, adapted from |Ouchi et al.H2010| l takes the Gunn-Peterson 
trough at z = 6.6 into account. For consistency with the analy¬ 
sis from lOuchi et al.| j2010^ , we use 2(7 limits in order to min¬ 
imise the number of potential interlopers (we show our 3(7 limits in 
Table 2). Finally, we use the NIR photometry to identify possible 
low-redshift interlopers which are extremely reddened by dust. For 
these sources, the Lyman-break is mimicked by the reddening from 
dust. However, it is possible to identify these interlopers based on 
their colours red-ward of the narrow-band. We exclude candidates 
which have J — K > 0.5 (empirically determined as a conservative 
lower limit) and are detected in the NIR by > 3a, which results in 


only 1 additional lower redshift interloper. This source, which is 
likely an [Oil] or [OIII]/H/I emitter, has an observed EW of 1050 
A. This additional step is not applied by Ouchi et al. i2010i, but 


it makes little difference for these luminosities. It is however im¬ 
portant for shallower narrow-band surveys, as will be shown in the 
following sections. The public SpUDS-IRAC catalogue (|Kim et al.| 
|2011[ | is based on a conservative magnitude limit and therefore does 
not contain faint enough objects, such that there is no match with 
any of the LAE candidates within a 2” radius, not even Himiko. It 
is however detected in deeper IRAC data ( [Ouchi et al.|20I3| ). 

We find a total of 99 Lyman-cr candidates in UDS and show the 
positions in the UDS panel in Fig.|^ The size of the symbols of the 
Lya candidates scales with the logarithm of the line-flux. 






























The Lya luminosity function at z = 6.6 1 


3.1.1 Lowering the EW criterion 

We retrieve 37 additional LAEs by lowering our excess (EW) cri¬ 
terion to 2 — NB921 > 0.5, but keeping the other conditions 
fixed. The risk of lowering the EW criterion is that the number 
of lower redshift interlopers increases. We use the NIR informa¬ 
tion and photometric redshift to exclude lower-redshift interlopers, 
because these will generally he very dusty (in order to mimic the 
Lyman-break) and thus bright and red in the NIR. In §3.1, we found 
only 1 interloper if we use an excess criterion of z — N31121 = 1. 
If we lower the criterion to z — NB921 > 0.5, we remarkably 
find only 4 interlopers, where we expected to find more. However, 
most of the additional 37 LAEs are very faint in the 2 -band, with 
magnitudes of ~ 25.5. This means that they need to have a very 
red colour (z — J > 1 or z — K > 2) to be detected in the 
near-infrared, since our detection limits are ~ 24 — 25 (see Ta¬ 
ble 2). Therefore, it is likely that a fraction of interlopers which do 
not have those extreme colours might be missed. Because of their 
low excess, the additional LAEs have faint luminosities and thus 
affect mostly the faint end slope in the luminosity function. How¬ 
ever, since we do not include these luminosities in the luminosity 
function due to their low completeness, the specific EW criterion 
used has little to no effect in our results. 


3.1.2 Comparison to Ouchi et al. (2010) 

The UKIDSS-UDS NB921 data has been analysed byjOuchi et al.| 
d^ot , who found 207 LAE candidates in UDS. The difference of 
108 in number of candidates arises since we are even more con¬ 
servative in our masked regions, limiting magnitude, and visually 
checking each candidate. However, when applying an analysis sim¬ 
ilar to Ouchi et al.’s, we find a very similar luminosity function (See 
§5.1 and §5.2). Erom the 16 spectroscopically confirmed LAEs by 
[Ouchi et al.| ( |2010l l, we have 15 LAEs in our first selection. The 
remaining source has a slightly lower excess in our analysis be¬ 
cause we do not use MAG-AUTO. Even though Lya is sometimes 
observed to be more extended than continuum emission (e.g. jStei-j 
[del et al.|[20TT]|Momose eraL|[20T4t , we choose to use aperture 
photometry for all our measurements. This is motivated by our line- 
flux completeness estimate in §5.1 and a comparison with detection 
completeness in §5.1.1. However, by lowering our excess criterion 
to 0.9, we also select this source. Since we recover the luminosity 
function and the spectroscopically confirmed LAEs, we find that 
our sample of LAEs is in agreement with the sample from |Ouchi| 
jet al.jpOTOt and that we fully recover their results. 

3.2 Selecting LAEs in COSMOS 

With the photometry presented in §2.4, we select line-emitters in 
the COSMOS field using: 

(i) z - NB921 > 1 

(ii) S > 3 

(iii) Pass visual inspection 

With the first two criteria, we find a total of 1633 candi¬ 
date line-emitters, although immediately exclude 1070 sources for 
which the excess is unphysically high (see Table 3 for the numbers 


Field 

Selection step 

Number of sources 

UDS 

Candidate excess sources 

1514 

0.81 deg2 

SpuriousAJn-physical 

147 

NB921 < 26.4 

Line-emitters 

1367 


Lyman-break selected 

100 


Red near-infrared 

1 


Final LAE candidates 

99 (15)* 

COSMOS 

Candidate excess sources 

1633 

0.96 deg2 

SpuriousAJn-physical 

1235 

NB921 < 25.8 

Line-emitters 

398 


Lyman-break selected 

19 


Red near-infrared 

2 


Variable sources 

2 


Final LAE candidates 

15 (2)* 

SA22-Deep 

Candidate excess sources 

359 

0.17 deg2 

SpuriousAJn-physical 

2 

NB921 < 26.4 

Line-emitters 

357 


Lyman-break selected 

6 


Red near-infrared 

4 


Final LAE candidates 

2 

SA22-Wide 

Candidate excess sources 

1674 

2.72 deg2 

SpuriousAJn-physical 

929 

NB921 < 24.3 

Line-emitters 

745 


Lyman-break selected 

25 


Red near-infrared 

5 


Variable sources 

1 


Final LAE candidates 

19 


Table 3. Number of sources in each selection step, for each field. The fi¬ 
nal numbers of LAE candidates are printed in bold. *The number between 
parentheses shows the number of spectroscopically confirmed LAEs. When 
comparing the fields, it can be seen that there are relatively more line- 
emitters selected as high-redshift source in UDS. This is because the UDS 
observations are deeper and the sample therefore exists of fainter observed 
sources, which generally are at higher redshifts. The number of interlopers 
identified based on near-infrared colours is relatively high in SA22, where 
the constraints on the Lyman-break are weaker due to shallower optical 
photometry. We could only check for variability in the COSMOS field and 
a small part of the SA22 field. 


after each step). We visually inspect the remaining ones and ex¬ 
clude a further 165. The number of candidates excluded based on 
an unphysical excess or visual checks is considerably higher than 
in UDS, which is caused by a high number of cosmic rays. The 
automatic rejection of cosmic rays is less successful because the 
number of raw exposures is roughly three times lower than in UDS. 
Other spurious sources have their flux boosted by haloes or spikes 
caused by bright stars. The difference with UDS is caused by the 
2 -band photometry, which is now relatively deeper compared to the 
narrow-band photometry. Our cut is slightly more conservative in 
this field than in UDS since we prefer completeness and robustness 
over the number of sources. We also rely on the UDS sample for the 
faintest luminosities (because of the deeper imaging) in our deter¬ 
mination of the luminosity function and use COSMOS for brighter 
sources. They agree in number densities (Eig.|^ indicating that this 
approach is correct. 

In total, we have 398 line-emitters. The photo-z distribution 
is peaked at [OHI]/H/3 at 2 = 0.83, with a smaller peak at Ha at 
2 = 0.40 and one around [OH] at 2 = 1.46. 10 have spectroscopic 
redshifts, of which 9 are [OHI]/H/3 and 1 [OH], see jSobral et al.j 
([2013^. In order to select LAEs, we apply the following criteria to 
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select high redshift line-emitters (see also §3.1): 

B > 27.9 A y > 27.3 ^ {i - z > 1.3 V i > 27.0) 

The optical limits are 2cr limits (also applied to UDS) computed by 
our empty aperture measurements, but consistent with |Muzzin et al.| 
( |2013^ . We use the deep near-infrared data to identify dusty lower- 
redshift interlopers with red near-infrared colours and identify two 
likely [OIII] or [Oil] emitters, which have J — K = 1.53 and 
J — K — 1.8 and observed EW of 420 A and 350 A, respectively. 

In addition to this, we are able to check our sources for 
variability. We have publicly available data from [Sobral et al.| 
( |2013^ which has been taken one year later than the COSMOS 
NB921 data and reaches a 3cr depth of ~ 25, which is similar 
to the magnitudes of our faintest Lya candidates. By comparing 
the magnitudes, we exclude the two brightest candidates (NB921 
~ 21.5), because they are completely undetected in the data from 
[Sobral et al.| ( |2013| l. This means that their brightness changes by 
~ 5 magnitudes and that they are thus likely supemovae. This 
means that similar surveys are expected to have roughly 2 SNe 
per 0.9 deg^ as contaminants to their sample of LAEs, except if 
the data has been split over multiple epochs. This may be very 
important to interpret results at higher redshift (e.g. [Faisst et al.j 
|2014[ [Matthee et al.|20f4) and shows how important it is to have 
data spread over time. 

We match our remaining LAE candidates to sources in the S- 
COSMOS-IRAC catalogue and find one match. This match is our 
brightest candidate, nicknamed “COSMOS REDSHIFT 7; CRT’. 
After these steps, we find 15 Lyman-o: candidates in COSMOS, of 
which two remarkably bright sources - brighter than Himiko - and 
show their positions in Fig. The size of the symbols of the Lyo? 
candidates scales with the logarithmic of the line-flux. 

We find no additional candidates when lowering the excess crite¬ 
rion. This is because the narrow-band is not as deep as in UDS, 
while especially the near-infrared constraints are stronger. There¬ 
fore, the exclusion of lower redshift interlopers is more successful, 
especially since we have shown in §3.1.1 that the additional can¬ 
didates are typically very faint, and therefore not present in our 
slightly shallower COSMOS imaging. 

3.2.1 Spectroscopic follow-up: bright sources in COSMOS 

In UDS, there are 16 spectroscopic confirmed LAEs by |Ouchi| 
[eTar] ( |2m0l l using Keck/DEIMOS. For COSMOS, we have ob¬ 
tained spectroscopic follow-up for our brightest two candidates 
using Keck/DEIMOS and DDT program 294.A-50I8 on VLT/X- 
SHOOTER and VLT/FORS2 presented in [Sobral et al.j j2015b[ l. 
Both of these are confirmed Lya emitters, at redshifts a = 6.604 
and 2 = 6.541, respectively. We nicknamed these galaxies CR7 
(see above) and MASOSZ0 CR7 is detected at only half of the 
narrow-band filter transmission, leading to a 2” luminosity of 
5.8 X 10"^® erg s“^, which is a factor 2 higher than our estimate 
from the photometry. Based on MAG-AUTO, the luminosity is 
9.6 X 10'*® erg s“*, and thus a factor 2.5 brighter than Himiko. 
This is a lower limit since the COSMOS NB921 observations are 
shallower than in UDS and might therefore miss some lower sur¬ 
face brightness regions and it assumes the «-band continuum be¬ 
ing flat. It is so extreme that is even detected individually in the 

® The nickname MASOSA consists of the initials of the first three authors 
of this paper. 


zYJHK bands and detected at 5 (j in the near-infrared stack, 
with YJHK = 24.9. It is also detected in IRAC, with a blue 
[3.6] — [4.5] colour, consistent with the [3.6] ^m flux being boosted 
by strong H/?/[OIII] line-emission (e.g. [Smit et al.|2014b[ l. Because 
of these broad-band detections, the source can be selected as a 
Lyman-break galaxy. Indeed, CR7 is present in the bright a ~ 7 
Ultra VISTA catalogs from [Bowler et al.j ( [20I2[ [2014^ . However, 
the Lya EW is much larger than the values used for the SED fitting 
(see [Sobral et al.||2015b[ i), meaning that the result from their SED 
fit requires revision and explaining the high x® • 

MASOSA (with a Lya luminosity of at least 3 x 10*® erg s“*, both 
in a 2” aperture and MAG-AUTO) is not detected in any of the 
near-infrared bands, and also not in the stacked image (meaning 
Y JHK > 26.7). MASOSA is undetected in the 2 -band, meaning 
that the luminosity is a lower limit. It is only brighter than Himiko 
when measured in 2” apertures, but this could also be due to our 
fainter narrow-band imaging in COSMOS than in UDS. It is not 
extended (diameter ~ 0.9”), while CR7 and Himiko show an ex¬ 
tent of ~ 3” in diameter. This means that the sources are of a dif¬ 
ferent nature and therefore interesting targets for follow-up study 
with e.g. HST. 

None of our LAE candidates are in the zCOSMOS jLilly et ^ 
[2009 [ l or UDSz ( [Bradshaw et al.[[2013[ [McLure et al.[[2013t cat- 
alogues, which we also did not expect, since they would likely be 
interlopers in that case. As expected, none of our LAE candidates 
except for CR7, are in the UV selected 2 ~ 7 galaxy catalogue 
from [Bowler et al.[ ( [2014^ , since they are not detected in broad¬ 
band photometry. None of the other sources in the [Bowler et al.[ 
( [2014^ catalogue are selected as line-emitters. 

3.3 Selecting LAEs in SA22 Wide & Deep 

Using our photometry described in §2.4, we select line-emitters 
using the following criteria: 

(i) z - NB921 > 1 

(ii) S > 3 

(iii) Pass visual inspection 

After these steps, we find 1674 line-emitters in SA22-Wide, 
from which we exclude 359 for having an unphysically high excess, 
and 347 emitters in SA22-Deep (see the corresponding panels in 
Fig. E and Table 3) and apply the following criteria to select high 
redshift line-emitters: 

u > 26.4 A g > 26.5 A r > 26.2 A (i — 2 > 1.3 V i > 25.9) 

The optical limits are 2a limits and the {i — z) criterion is similar 
to that in the other fields. We visually check each LAE candidate 
in the NB921 images and exclude a further 560 in SA22-Wide 
and 2 in SA22-Deep. Most of these objects are cosmic rays which 
were not detected automatically and sources which have their flux 
boosted by haloes or spikes caused by bright stars and are thus 
not real excess sources. In SA22-Wide our E-value from pixel 
based measurements corresponds to a 3 E if we measured the 
background with empty apertures, and to 1.5 E in SA22-Deep. 
These differences are caused by a different limiting narrow-band 
magnitude. 
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Figure 3. Position on the sky of the three survey fields. The three panels show the relative areas of the three fields. Grey dots show all detections with 
A^S921 < 22 (chosen in order to control the file size of the images), highlighting the masked regions due to bright stars as empty circles and masked noisy 
areas which are due to our pointing strategy (such as the grid pattern in SA22). In SA22, the detections in the Deep region have a slightly darker colour. The 
green symbols show the positions of our LAEs. The size of these symbols scales with the line-flux (luminosity), following size ~ (logio(LLyQ,))^/^, and has 
the same limits and range for all the three fields. From this, it can be seen that UDS is significantly deeper since it has many more LAEs with small symbols. 
We furthermore note that two LAE candidates in SA22 are not visible in the image because they overlap with other symbols due to a small separation on the 
sky (~ 4”). 
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There is one LAE candidate in SA22-Wide which is in the small 
overlapping region with SA22-Deep. This source happens to be 
variable, as it is undetected in the SA22-Deep data. We could not 
check the major part of the SA22-Wide field for variability and 
therefore use a statistical correction to the luminosity function 
using our empirical results from the COSMOS region, where we 
found 2 variable sources (likely supemovae) per 0.9 deg^. The area 
is 2.72 deg^, so we weight our number densities down by 6 sources. 

Since the optical photometry in SA22 is shallower than in the 
other fields (see Table 2), there is a higher chance of our candidates 
being lower redshift interlopers, since the Lyman-break constraints 
are not as stringent. Using the near-infrared J and K data, we 
identify objects with significantly red NIR colours (J — K > 0.5). 
In SA22-Deep, we exclude 4 out of our 6 LAE candidates based 
on optical colours only, and are thus left with 2 LAEs. These four 
lower redshift interlopers can be called extreme emission line 
galaxies, since their observed EW is ~ 400 A. In SA22-Wide, 
we exclude 5 interlopers out of the 24 LAE candidates from the 
near-infrared photometry. Because the SA22-Wide candidates are 
brighter, possible interlopers are easier to exclude based on their 
optical colours. Therefore, this number is relatively lower than in 
SA22-Deep. In total, we find 2 LAEs in SA22-Deep, and 19 LAE 
candidates in SA22-Wide. Their spatial position is shown in Eig. 

where the size of the symbols scales with the logarithm of the 
line-flux. 

The LAE candidates in SA22-Wide are particularly bright, 
with luminosities 3 — 16 x 10^® erg s“^, if they are at z = 6.6. 
We note that 4 candidates are in pairs which are separated only by 
~ 3 — 5” in the sky (such that only one point is seen in the Eig.|^. 
Once these candidates are spectroscopically confirmed, this sample 
will allow us to study the variety of bright LAEs. 


4 NUMBER COUNTS, COMPLETENESS AND 
CORRECTIONS FOR FILTER PROFILE BIAS 

Our main diagnostic is the luminosity function and its evolution 
with redshift. The luminosity function shows the volume number 
density of Lya emitters with a certain (observed) luminosity. In 
this section, we derive the observed number counts and all other in¬ 
put values for the luminosity function. This is because in addition 
to the raw number counts per bin, there are important corrections 
to be made, since our observations and analysis introduce biases 
and systematic errors, for which we correct in the following sub¬ 
sections. 

The probed volume can be calculated using the FWHM of the 
filter, since that puts a lower limit and an upper limit to the Lya 
emission-line redshift, which are 6.50 and 6.61 respectively. We 
calculate the volume then as the difference in comoving spherical 
volumes within the upper and lower redshift limits, and multiply 
this by the fraction of the sky that is our survey area. We find a 
comoving volume of 9.02 x 10® Mpc® deg“^, corresponding to 
7.4 X 10® Mpc® in UDS, 8.7 x 10® Mpc® in COSMOS, 1.5 x 10® 
Mpc® in SA22-Deep and 2.5 x 10® Mpc® in SA22-Wide. In total 
our volume is 42.6 x 10® Mpc®. 


Field 

80% completeness flux 

UDS 

1.3 X 10“^^ ergs“^ cm“^ 

COSMOS 

4.4 X 10“^^ erg cm“^ 

SA22-Deep 

7.4 X 10“^^ erg cm~^ 

SA22-Wide 

10.0 X 10“^^ erg cm“^ 


Table 4. The line-flux for which our completeness is 80 %, shown in our 
different fields. This depends on both the NB921 and z band depths. Note 
that, because of this, even though the narrow-band photometry has similar 
depths in SA22-Deep and UDS, the line-flux completeness is very different 
due to a different broad-band limit. This highlights the need for line-flux 
completeness over detection completeness. 

4.1 Line-flux completeness 

Our selection of line-emitters relies on the measured narrow-band 
excess and the excess significance. This means that photometric 
errors can lead to missing real Lya emitters at z = 6.6, especially 
at the faintest luminosities. The result is that our sample is 
incomplete. How incomplete our survey is at a given luminosity 
(line-flux) depends on the survey depth, source extraction and 
selection method. We can measure the incompleteness with a 
simulation based on a sample of observed sources which are 
consistent with being at high redshift (z > 3, using Lyman-break 
criteria), but are not detected as line-emitters. This sample is 
detected and analysed in exactly the same way as our sample of 
Lya candidates and has similar narrow-band magnitude distribu¬ 
tion. To these sources (> 1,000 in total in the three deep fields, 
and > 10,000 in SA22-Wide) we artificially add line-flux by 
changing their NB921 and z-band magnitude correspondingly, 
and test whether it is then selected as a line-emitter based on 
the updated excess and excess significance. The completeness 
is then obtained for each line-flux by measuring the fraction of 
sources which is selected as line-emitter after adding the flux (e.g. 
ISobral et al.||201^ . We measure the completeness in the three 
deep fields for each pointing, but averaged over the SA22-Wide 
pointings. Because of the limited depth of SA22-Wide, there are 
less number of sources per pointing and the statistics are thus 
weaker. We find that the completeness in different pointings in 
UDS and COSMOS are very similar. We show the line-flux for 
which each of our fields is complete up to 80 % in Table 4. 
These fluxes correspond to luminosities from 6.37 x lO'^^ erg 
s”^ (UDS) up to 4.9 X 10^® erg s"® (SA22-Wide). Note that 
even though the narrow-band photometry has similar depths in 
SA22-Deep and UDS, the line-flux completeness is very different 
due to a different broad-band limit. This means that a measure of 
completeness based on detection only will give inconsistent results. 

For each luminosity bin, we correct the number of sources 
by dividing by the completeness. We also divide the poissonian 
errors by this completeness value. Only bins with a completeness 
higher than 40 % are included in the luminosity function. The 
completeness correction is strongest for the faintest luminosities, 
and thus increases the number density mostly at low luminosities. 


4.1.1 Detection completeness in SA22 

In order to access the quality of our data in our wide SA22 survey, 
we use our three spectroscopically confirmed bright LAEs to 
estimate the detection completeness. By placing them at random 
positions in our images and see whether we recover them, we 
know whether our data is sufficient to observe these luminous 
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Figure 4. Number counts in our three fields, compared to the bins from 
lOuchi et al.H2010y The bins are only coirected for completeness. Our bins 
in UDS vary with this from |Ouchi et al.|j2010^ because we use luminosities 
based on 2” apertures and apply a different completeness correction. The 
SA22-Wide bins are corrected for contamination from variable sources and 
supernovae, empirically calibrated in COSMOS. Note the good agreement 
between measurements from all fields, although some variance is found 
which is likely due to cosmic variance. 

sources and we can furthermore check our line-flux completeness 
procedure (see above) and see how the two compare. 

We produce small cutouts (5” x 5”) around CR7, Himiko and 
MASOSA and add them to 100 random positions per pointing in 
SA22, excluding masked regions. After this, we run SEXTRACTOR 
with identical settings as used on the original images and compute 
the fraction of our input sources which is detected. We repeat this 
a 1000 times per image and use the average recovered fraction 
as detection completeness. On average, we find a detection com¬ 
pleteness of 44 %, with a standard deviation of 20 % in different 
pointings. The detection completeness is highest for MASOSA, 64 
%, around the average for CR7, 43 %, and lowest for Himiko, 27 
%. This is because the first source is not extended, while the other 
two are extended and therefore have lower surface brightnesses. 
Note that we do not exclude pixel positions with actual sources 
or regions with a slightly lower signal to noise, which both 
decreases the completeness. The average detection completeness 
is remarkably similar to our estimated line-flux completeness 
(which is 46 % for the average line-flux of the three sources). The 
large variation in detection completeness between the different 
sources, which have almost the same 2” magnitude, highlights 
the need for a completeness based on line-flux, instead of detection. 



logio Uya (erg 

Figure 5. Luminosity function at z = 6.6: comparison of number densities 
and fits with and without filter profile correction. We compai'e our number 
densities in UDS with those of |Ouchi et al.| i [2QTQ) , which are largely based 
on the UDS field as well. We show that our red hexagons (before correcting 
for the filter profile) agree well with the green squares from |Ouchi et al.| 
jlQlO) , whose fit to the data is shown as a solid green line. The dashed 
green line shows our fit to their total data (fixing a = —1.5 in Eq. 5), 
which differs significantly from their published fit. The fit to our data (a 
fixed to —1.5; dashed red line) again agrees well, indicating that our results 
are similai*. The effect of the filter profile con'ection is shown by compar¬ 
ing the blue hexagons with the red hexagons. The effect is that the number 
density of bright line-emitters is higher, while the number density of faint 
line-emitter is slightly lower. The blue line shows the fit to the bins after cor¬ 
recting for the filter profile, which again highlights the effect of the correc¬ 
tion. The grey line shows a model prediction by |Gronke et al.|j2015) which 
is based on the LEG LF and a Lyo EW distribution, frozen at z = 6.0. It is 
remai'kable that there it agrees well with the blue curve, despite not being a 
fit. 


Luminosity bin 

Number density correction factor 

42.5 

0.99 

42.7 

1.07 

42.9 

1.18 

43.1 

1.32 

43.3 

1.51 

43.5 

1.77 

43.7 

2.08 

43.9 

2.79 


Table 5. Correction factors for the number densities at z = 6.6. These 
corrections are made for the bias arising from the observations through the 
filter profile not being a top-hat. Because of the filter profile, luminous LAEs 
can be observed as faint LAEs, meaning that their real number densities are 
higher than observed. This is particularly important for when comparing 
narrow-band LAE searches with IFU based LAE searches. 


4.2 Number densities 

We show our number densities in Fig. [Sand compare with the num¬ 
ber densities from |Ouchi et al.| l |201^ (purple circles), which is 
based majorly on UDS. Our UDS points agree with those of|Ouchi| 
|et al.| l [20lQl t, while the SA22-Deep and COSMOS bins (which are 
spectroscopically confirmed) converge at brighter luminosities and 
are also consistent with |Ouchi et al.| ( [2010| t. Our SA22-Wide num¬ 
ber densities are more uncertain, since there is no spectroscopic 
confirmation yet and the photometric constraints are weaker than 
in the other fields. However, even if there are still some contami¬ 


nants, these further highlight a departure from a Schechter function 
(already indicated by our spectroscopically confirmed sample) at 
high luminosities and indicate that the observed Lya luminosity 
function at z = 6.6 can be fitted by a powerlaw (e.g. the pentagons 
in Fig.^. The powerlaw fit is: 


logio ( 


-I- 

Mpc“^ 


) 


68.38 - 1.681ogio(^^l^) 
ergs“^ 


Since we have only two sources in SA22-Deep and since this 
agrees very well UDS and COSMOS, we will include them when 
we refer to the UDS-l-COSMOS sample in the remainder of the 
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Data — set 

UDS + COSMOS + SA22 

UDS + COSMOS 

UDS 

UDS + COSMOS + SA22 without filter correction 

UDS without filter correction 

Ouchi et al. 1^20101 data, our fit 

UDS + COSMOS + SA22 

UDS + COSMOS 

UDS 


UDS + COSMOS + SA22 without filter connection 
UDS without filter correction 
Ouchi et al.|i 20I0| data, our fit 


g [fixed] logio(<I>*)[Mpc logio(L*)[ergs 


-2.0 

-4 

^■'^'^-0.12 

43 61+^-^^ 

06 

3.14 

-2.0 

-4 13+0-10 

43 31+0-09 

2.49 

-2.0 

-4 

^■^°-0.44 

43 34+0.38 

1.78 

-2.0 

-4 40+0-10 

43 42+0.10 

1.75 

-2.0 

.3 0*7+0.15 

H-j.iZ-o 09 

1.31 

-2.0 

0 70+0.13 
-i./5 o 18 

^2».uo o 08 

0.48 


-1.5 

0 QT+U.Ub 

-•J-^-J_0.05 

43 33'''^'^'^ 

6.65 

-1.5 

.3 ^ 2 + 0-06 

43.05t°;°® 

2.07 

-1.5 



0.95 

-1.5 

.3 01 +0.06 

43.20t°;°i 

3.87 

-1.5 

o c*5+0-09 
-•J-J-J_0.10 

42.88t°:°^ 

0.68 

-1.5 

O '7C+0-08 
-■:>-JJ 0.08 

42.84^0;°^ 

0.71 


Table 6. Values to our Schechter fits to the luminosity functions. Because of our limited depth, we fix the faint end slope to either —2 or —1.5 and only keep 
<&* and L* as free parameters. 


Data-set 



a 



b 

xL 

UDS -H COSMOS - 1 - SA22 

62.07 

± 

3.45 

-1.53 

± 

0.08 

2.75 

UDS -H COSMOS 

60.84 

± 

5.00 

-1.50 

± 

0.12 

4.72 

UDS + COSMOS + SA22 before filter correction 

68.38 

± 

4.58 

-1.68 

± 

0.11 

1.62 

UDS + COSMOS before filter correction 

66.26 

± 

6.10 

-1.63 

± 

0.14 

2.67 


Table 7. Power-law fits to the number densities of the functional form: logj^Q () = a -|- b log jg () The first two are a fit to the corrected number 
densities (Fig.[^ and the second to the observed (Fig.|^ number densities. We show both fits with including all three fields, and only for UDS and COSMOS. 


text. We will also refer to the SA22-Wide results as SA22. 

Since our LF estimate is based on binning the data, we suffer 
from Eddington bias (c.f. |Ilbert et al.|2013) . As luminosities have 
photometric errors, these uncertainties scatter sources from one 
bin to the next bin. However, due to the shape of the LF (more 
sources in the fainter bins than in the bright bins), the luminosity 
uncertainties move more sources into the luminous end than vice 
versa. Therefore, this bias tends to overestimate the number of 
luminous sources and underestimates the fainter sources. It could 
therefore overestimate the bright end. We are aware of this bias, 
but do not apply a correction in order for consistency with previous 
surveys and since its effect are similar between different redshifts. 


4.3 Filter profile bias correction 

Since our filter is not a perfect top-hat, the exact redshift of the Lyo? 
emission line influences the observed luminosity. This means that 
intrinsic luminous LAEs which are detected at the edges of the fil¬ 
ter (where the transmission is lower) are observed as fainter LAEs. 
It also means that the probed volume depends on the luminosity, 
since luminous sources can be detected over a larger redshift 
range, but will be observed as fainter sources. For example, our 
brightest spectroscopically confirmed source in COSMOS, CR7, is 
actually detected at only 50 % of the transmission and is thus even 
brighter than our photometric estimate. Corrections for this effect 
are derived with a simulation, similar to |Sobral et al.| ( |2013| l for 
Ha line-emitters. We use the Schechter fit of our UDS-l-COSMOS 
data to generate a million Lya emitters and assume that they have 
a random redshift between the edges of the filter. We then convolve 
the luminosities with the filter profile into an observed population. 
Corrections are then obtained by comparing the number of sources 
in each luminosity bin before and after applying the filter-profile. 
The result of our correction is that the number-density of luminous 
sources is increased, while it decreases at low luminosities. We 


note explicitly that this correction is required to remove the bias 
from observation strategy, since e.g. an IFU survey without a filter 
would not suffer from this bias, and that it is not related to any 
intrinsic effect of the sources. We show the correction factors for 
2 = 6.6 in Table 5. 


4.3.1 Comparison to \Ouchi et al.\^20W^ and effect of the filter 
profile correction 

There are two reasons for the small differences between our 
UDS and Ouchi’s number densities (Fig. |^: the first is that our 
completeness correction is based on line-flux and the selection 
of emitters, while their completeness correction is based on 
detection completeness and narrow-band magnitudes. The other 
difference is that [Ouchi et al.| P010| l uses MAG-AUTO to 
estimate narrow-band magnitudes (which are used to compute 
luminosities), while we use the magnitude in 2” apertures. As 
Lya is often extended, MAG-AUTO might have been a better 
choice. We however chose to consistently use the 2” aperture 
since we also used this for our (important) line-flux completeness 
correction. Using MAG-AUTO would mean the introduction of 
an additional uncertainty. We compare the luminosities derived 
with both 2” apertures and MAG-AUTO for the spectroscopic 
confirmed sample in UDS and find that corrections of -l-0.11 dex 
can be used to statistically correct the luminosities. This is used in 
Fig. |4] and in all other following luminosity functions. Our results 
do not strongly depend on this correction. 

The effect of the filter profile correction on the luminosity 
function derived by jOuchi et al.| ( [2010[ l is shown in Fig.[^ We first 
fit a Schechter function to our UDS data, both before and after cor¬ 
recting for the profile (respectively the dashed red and solid blue 
line, see also Table 6) and then compare this fit to the Schechter 
function from [Ouchi et al.j p010[ >. Because the brightest bin in 
their sample contains only one source (Himiko), the error on the 
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Figure 6. The Lyman-a luminosity function at 2 = 6.6 compared to literature data. Our most robust luminosity function is shown as a solid blue line. This is a 
Schechter fit to our combined UDS and COSMOS data (blue circles, see also Table 6), for which the brightest LAEs have all been confirmed spectroscopically. 
Our additional SA22 data is shown in open circles and is consistent with the upper limits from our robust sample. We also place upper limits (blue and open 
arrow) at the luminosity bin just brighter than the most luminous obseiwed sources, meaning that there is less than one of these in the probed volume. The 
dashed blue line is our power-law fit (see Table 7) to the data from all three fields. The fit from |Ouchi et al.|(2QlQ) at z = 6.6 differs for two main reasons (see 
also Fig. [ 5 ), namely practically not including the brightest bin to their fit (due to very large errors, as the fit contains only a single source) and not correcting 
for different biases caused by the filter profile. This is also the major reason why our results are slightly different with the results from|Kashikawa et al.|p006| 
(ochre diamonds) and |Hu et al.H201^ (red squares). Other reasons are cosmic variance, since they only probed small areas (Kashikawa et al.), and small 
(spectroscopic) completeness (Hu et al.). 


bin is extreme. Therefore, the fit from Ouchi et al. predicts a 30 
times lower number density of bright LAEs than observed. We also 
fit a Schechter function to their data in log space (green dashed 
line). Note that we fix a to —1.5 in all our fits, similar to previ¬ 
ous searches (e.g. [Ouchi et al.|2008] l, since even the deepest UDS 
data is too shallow to constrain a accurately. However, in Table 6 
we also provide fit values for an a fixed to —2. From Fig. it is 
clear that, first of all, our luminosity function is in agreement with 
[Ouchi et al.| ( |2010l l if we use similar corrections and include their 
brightest bin to the fit (as can be seen by comparing the two dashed 
lines in Fig.[^. Second, the effect of the filter profile is highlighted 
by comparing the solid blue to the dashed red line; after correcting 
for observational biases from the filter profile, the number density 
increases mostly at brightest luminosities. 


5 LYa LUMINOSITY FUNCTION AT Z = 6.6 

In this section, we present the 2 = 6.6 Lya LF from our combined 
analysis in UDS, COSMOS and SA22. As a functional form, we 


use the well-known Schechter function: 

(j){L)dL = 4>*{L/L*Y exp(-L/L*)d(L/L*) (5) 

We convert our observed line-fluxes to luminosities by assuming a 
luminosity distance corresponding to a redshift of 6.56, which is 
the redshift of the center of the filter. We combine the luminosities 
in bins with widths of 0.2 dex and count the number of sources 
within each bin and correct this number for incompleteness. The 
errors on the bins are taken to be Poissonian. The number of 
sources is divided by the probed volume, such that we obtain 
a number density. We then apply our corrections for the filter 
profile bias. Only data where the completeness is at least 40 % 
is included. The resulting luminosity function is shown in Fig. 
where we also compare with other published 2 = 6.6 LAE data. 
The evolution between 2 = 5.7 and 2 = 6.6 is shown in Fig.[^ 
while the left panel of Fig.j^shows the evolution towards 2 = 7.3. 
We are cautious about interpreting the results from SA22 because 
of the less stringent photometric criteria, even though they fully 
agree with results from the other fields. The results from UDS and 
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COSMOS however, are confirmed by spectroscopy. 

It is interesting to compare these results with the model from 
IGronke et al.| ( |2015^ , which is shown as the black line in Fig. 
This model uses the UV LF and a probability distribution (PDF) 
of Lyof EWs to predict the Lya LF. The EW distribution generally 
evolves with redshift, but in this case, it is frozen to the EW 
PDF at 2 = 6.0, because of possible effects from re-ionisation. 
It is remarkable that the prediction from |Gronke et al.| j2015^ 
seems to be consistent with our blue points. Differences arise 
because of their steeper faint end slope (~ —2.2), which is largely 
unconstrained by the depth of our current data of LAEs. The 
agreement highlights the need for the correction of the filter profile 
bias when comparing narrow-band derived LFs with LFs derived 
from spectroscopy (either follow-up or blind IFU). 


As noted in §5.2.1, our results in UDS differ by those from 
|Ouchi et al.| j2010^ at brighter luminosities due to a different 
treatment of the brightest bin in the fit (solid green line) and by 
correcting for the filter profile (which effect is shown in Fig. |^. 
This explains also the differences (although largely within the 
errors) with |Kashikawa et al.| ( |2006^ , although cosmic variance 
plays a role because of their limited survey area. As noted by 
Kashikawa et al.| (|201 l[l, the difference between the results from 


Hu et al. (|2010|) and the others is due to incompleteness of the 


sample of Hu et aL| ( |2010^ , since they rely on spectroscopic 
follow-up with too short integration times. Before correcting 
for the filter profile, our results in COSMOS agree with those 
from [Ouchi et al.| ( |2010| l (see Fig. |^, but even after correcting 
for the filter profile, our combined UDS-l-COSMOS LF agrees 
with the brightest bin of |Ouchi et al.|(2010| >, and disagrees only 
slightly with the second brightest bin. Our SA22 results are not 
yet confirmed spectroscopically and are thus upper limits when 
viewed most conservatively. There is however excellent agreement 
with the spectroscopically confirmed COSMOS sources and with 
the upper limits from the UDS-l-COSMOS sample (blue arrow in 
Fig.|§l. If all (or even only a fraction) of these very bright Lyo? 
emitters in SA22 are confirmed, this indicates that the observed 
Lya luminosity function at z = 6.6 can be fitted by a powerlaw 
(e.g. the pentagon s in Fig.[^, similar t o the UV luminosity function 
at 2 = 6 — 7 (e.g. |Bowler et al.|2014| l. 


Now we have used a combination of wide and ultra-deep 
fields, we have established a new LF at 2 = 6.6 (see Table 6) and 
will compare this with results in the literature at lower redshift. To 
be conservative, we will use our fit to the UDS-l-COSMOS sample 
for comparison, although our results are only strengthened by the 
results from SA22. 


5.1 The z = 5.7 — 6.6 evolution of the LF 


Previous studies (e.g. van Breukelen et al.|2005 Ajiki et al. 2006| 


IShimasaku et al.||2006| [Gronwall et al.||2007[ [Ouchi et al.| 20081 

have found that the observed Lya luminosity function is remark¬ 
ably constant between 2 = 3 — 6. While the LBG LF declines over 
this redshift range, the implication is that the strength of the Lyo? 
line increases with redshift, which is also confirmed spectroscopi¬ 
cally (e.g. |Stark et al.|20T0l l. At 2 > 6, there is evidence both from a 
declining success-rate of spectroscopic observations of LBGs (e.g. 
Schenker et al.|2012[l, and also for a drop in the Lya LF l |Ouchi| 


et al.poiof Konno et al.||2014^. In our independent analysis we 


confirm evolution of the Lya LF from 2 = 5.7 to 2 = 6.6, al¬ 


though only robustly at fainter luminosities, L < 10'*® erg s *, see 
Fig.[7| As a comparison at 2 = 5.7, we use data-points from |OucTu] 
|et al^2008t (green squares in Fig.l^, complemented with the wide 
area survey from| Westra et al.| l| 20()^ . We were not able to compare 
our results to the number densities of |Murayama et al^ l |2007| l at 
2 = 5.7 (2 deg® in COSMOS) since they do not correct for com¬ 
pleteness. As discussed above, we apply a correction for the filter 
profile, based on a generated sample of a million LAEs following 
the 2 = 5.7 Schechter function from |Ouchi et al.|j2008t . The result 
is shown as the orange squares in Eig.[^ For the data from|Westra| 
|et al.| l (^06| l, we only show the corrected bins. We then fit a new 
Schechter function to the corrected points, shown as a solid orange 
line. The evolution between 2 = 5.7 and 2 = 6.6 can be seen by 
comparing the orange and blue solid lines. It can be seen that at 
fainter luminosities (hLya ~ 10*®'® erg s“*) the observed number 
density of LAEs declines by a factor ~ 3, while there is no evolu¬ 
tion at the bright end hLya > 10*®'® erg s“*). Note that we only 
compare to our 2 = 6.6 sample based on UDS and COSMOS and 
that the results at the bright end are confirmed spectroscopically. 
Addition of the SA22 bins only strengthens the conclusion. Unfor¬ 
tunately, however, it is not yet possible to compare the SA22 results 
with lower redshifts, since there is no comparable data-set at 2 < 6. 
Even without the filter profile correction, we find differential evo¬ 
lution of the luminosity function. We will discuss explanations for 
this differential evolution in the discussion in §6.1 and §6.2. When 
adding even higher redshift data from 2 = 7.0 jlye et al.|2006 | |0t^ 
|et al.|20l0) and 2 = 7.3 jShibuya et al.|20l2 JKonno et al.|2014 1 
narrow-band surveys (see the left panel of Eig. y, the LF keeps de- 
creasing at luminosities smaller than 10*® erg s“*. Note that these 
bins consist of just a handful of sources and that we have not cor¬ 
rected them for the bias due to the filter profile. As the survey areas 
of these surveys are below 0.3 deg®, they are all limited severely 
by cosmic variance and it is therefore not possible to compare the 
evolution at the bright end. A few additional shallower pointings 
will be useful to constrain the evolution of the bright end as well, 
and confirm or neglect a continuing differential evolution. 


6 DISCUSSION 

6.1 LF evolution and re-ionisation 

The evolution of the observed Lya LF can be caused by different 
processes: 

i) an intrinsic dimming of the Lya EW. As the cosmic star forma¬ 
tion rate density derived from UV observations at 2 > 6 continues 
to decline with redshift (e.g. |Bouwens et al.|2014b[|McLeod et al.| 
|2014| l, the production rate of Lya can decrease due to a lower star 
formation rate, since Lya is emitted by the recombination of hydro¬ 
gen atoms which are photo-ionised by massive, young, short-lived 
stars. However, the cosmic star formation rate is already declining 
from at least 2 > 3, while the observed Lya LE is constant be¬ 
tween 2 = 3 — 6, indicating an increase in Lya EW with increasing 
redshift (see also |Stark et al.lpOlO) ) and it is hard to explain a sud¬ 
den reversal of this trend. Eurthermore, the star formation towards 
higher redshifts can be partly contributed by formation of metal 
free Pop III stars, which are predicted to produce copious amounts 
of Lya emission l |Schaerer|2003^ . 

ii) the observed drop in the Lya LF can be explained by a 
lower escape of Lya in the interstellar medium of galaxies l|Dayal| 
|& Ferrara|2012[[Dijkstra et al.|2014^ . This escape fraction is largely 
unconstrained, although first direct measurements of matched Ha- 
Lya observations at 2 = 2.2 indicate an average escape of 5 % 
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Figure 7. Evolution of the Lya luminosity function from = 6.6 to 2 = 5.7. We compare our 2 = 6.6 LF (blue solid line) to published data and fits to the 
data at 2 = 5.7. We apply a filter profile comection to the data from|Ouchi et al.|J2008| at 2 = 5.7 (green: uncorrected, ochre: corrected) and fit a Schechter 
function to the corrected data (orange solid line). We also show the spectroscopically confirmed results from |Westra eljilj j20()6|at 2 = 5.7. Comparing the 
orange solid line to the blue solid line shows that the LF evolves only at the faint end and not (as has been claimed by |Ouchi et al.lpOlO) ) at all luminosities. 
In fact, we find no evolution for L > 10^®'® erg s"'^. We show in §6.2 that this may be a consequence of re-ionisation. We do not show our fit which includes 
the SA22 data-points, since there is no comparison available at 2 = 5.7. This motivates the need for larger volumes at 2 = 5.7 as well, in order to see if there 
is evolution at the bright end. 


([ Hayes et al.|2010| ). However, other, more indirect measurements of 
the escape fraction in LAEs is in general higher (30 %, e.g.jWard-j 
[low et al.|2014[[Kusakabe et al.|2014) . Using a joint analysis of the 
evolution of the Lya and LEG luminosity functions, [Hayes et al.j 
( [201 1^ finds that the volumetric (statistical) Lya escape fraction 
increases with redshift up to 2 = 6 and decreases at higher red- 
shift, although this result is based on the integration of the Lya LF 
from [Ouchi et al.[ ( [2010[ (. As we show, in this paper, their large er¬ 
ror on the bright end results in an underestimate of the luminosity 
density. The intrinsic Lya escape fraction can be measured directly 
at 2 > 3 once infrared spectroscopy of the Ha line with the James 
Webb Space Telescope is possible. It is thought that the primary 
factor in driving the escape fraction is the neutral hydrogen column 
density (e.g. [Nakajima et al.[|2012[ (, but kinematics and dust pro¬ 
duction can also have a role. Similar to the intrinsic dimming of 
Lya EW, it is hard to explain a reversal of the trend highlighted 
by [Hayes et al.[ ( (201 1[ |, i.e. an increase up to 2 ~ 6 and a decline 
afterwards. 

iii) a sudden decrease in the observability of Lya can possi¬ 
bly be explained by an increase in the incidence of dense pockets 
of neutral hydrogen in the line of sight (e.g. Lyman Limit Systems 


or Damped Lyman-a absorbers). [Bolton & Haehnelt[j2013| l argue 
that the majority of optically thick gas along lines of sight is present 
in absorption systems, although again, there is no definite explana¬ 
tion why there should be a sudden increase, or why this increase 
depends strongly on Lya luminosity. 


iv) the most quoted reason (e.g. Santos|2004[ Haiman & Cen 

2005 

|McQuinn et aL||2007| |Choudhury et al.|2014| Jensen et al. 

2014 

1 is that the observed decrease in both the evolution in the Lya 


LF and the declining spectroscopic success-rate in LEG follow¬ 
up, is caused by an increased neutral fraction of the IGM, which 


leads to a higher opacity to Lya photons. [Mesinger et al.[ j2015[ l 
however argue that re-ionisation can not be the only driver of the 
observed changing opacity to Lya and [Konno et al.[ ( |2014^ mention 
that the neutral IGM fraction required to explain their observed 
evolution of the Lya LF (see left panel of Fig. is in tension 
with results from the polarisation of the Cosmic Microwave Back¬ 
ground (CMB) ( [Planck Collaboration et al.|2014^ , with a mean re¬ 
ionisation redshift of 2 = 11.1. On the other hand, [Choudhury et al.[ 
( [2014^ argue that the most recent CMB results are not in disfavour 
with a re-ionisation epoch that ends between 2 = 6 — 7. The latest 
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Figure 8. Left: Evolution of the luminosity function evolution from z = 7.3 to 2 = 5.7. We compare our z = 6.6 LF (blue solid line) to published data 
at 2 > 7|Konno et al.||2014| (red squares)|Shibuya et al.||2012|;|Ota et al.|p010);|Iye et al.|^2006J (purple triangles, pentagons and diamond, respectively). 
We also show our LF fit to the corrected 2 = 5.7 data (orange solid line) and the green squares show the number densities at 2 = 5.7 from|Ouchi et al.| 
j2008) . Note that the errors on the 2 > 7 data are still significant and that the surveys are limited to small areas, meaning that the LF is still unconstrained at 
luminosities > 10“^®'^ erg/s. A few additional, shallower pointings with the Subaim S-Cam and the NBlOl (z = 7.3) filter would place useful constraints on 
the evolution of the Lya LF at these epochs. Right: Toy-model evolution of the luminosity function in a neutral IGM. The black line shows the input LF, which 
is the 2 = 5.7 LF from |Ouchi et al.|j2008) . We fix the parameters controlling the age and number of escaping ionising photons from LAEs (to age =100 Myr, 
fesc,ion = 5% and fesc,Lya = 30 %) and investigate how a changing neutral fraction (Xuj) influences the observed LF. As can be seen, the evolution starts 
at the faintest luminosities only and gradually towards higher luminosities. The highest luminosities (> 10“^® ® erg s“^) are still observable. Faint LAEs are 
observable if they are in the ionised sphere of a neighbouring, luminous source. Therefore, the observed faint-end slope depends on the clustering. The turnover 
luminosity, at which LAEs can ionise their own surrounding enough, depends on our input parameters and the neutral fraction. Note that our specific values for 
the neutral fraction are arbitrary since our input parameters are largely undetermined, and they are only showed for illustration purposes. While typical models 
require much larger changes in the neutral fraction (e.g. |Dijkstra|^2014) ), our values are particularly low because we assume that all Lyo emission is absorbed 
when the ionised sphere is < 1 Mpc, while in reality only part of the line-emission is absorbed, and we ignore other processes such as galaxy outflows which 
can increase the observability of Lya. 


Planck results suggest a re-ionisation redshift of 2 = 8.8 ( |Planck | 
ICollaboration et al.|2015| >, but this result is model dependent. 

A lot of these explanations are very degenerate with the cur¬ 
rently available observational data and can only be solved by disen¬ 
tangling the effect from galaxy formation (e.g. a varying intrinsic 
production or escape of Lya photons) with a cosmological cause 
(re-ionisation of the IGM). One way to overcome this is to study 
the clustering properties of LAEs in wide area surveys, as is the 
goal of the HyperSuprimeCam survey on the Subaru telescope (see 
also §6.2). [jTOsen et al^ ( |2014^ show that an increased (observed) 
clustering signal of LAEs is indicative of a more neutral IGM, be¬ 
cause Lya is preferentially observed in an ionised region and more 
clustered LAEs will lead to a larger ionised region around them. 
Mimicking this clustering signal is unlikely by changes in the in¬ 
trinsic escape fraction or production of Lya. Another measurement 
which will overcome the degeneracy is to measure the intrinsic Lya 
escape fraction directly using spectroscopic Ha measurements, for 
example possible with a matched narrow-band survey to Lya with 
a wide field camera which is sensitive to 2-5 micron radiation in 
space (since the night sky background is very high in these wave¬ 
length regions). 

In the next section, we argue that our observed differential 
evolution of the Lya LE can be explained by a simple model for 
re-ionisation, although we do not claim that this is the only possi¬ 
bility. 


6.2 Toy-model for the evolution of the Lya LF during 
re-ionisation 

Since the process of re-ionisation likely has a patchy nature, where 
the Universe is ionised either bottom-up (overdense, early collaps¬ 
ing regions first) or top-down (voids first), there is an observable 
effect in the evolution of the Lya LF. |Haiman & Cen|p005| > argue 
that the observability of bright LAEs is less attenuated by a neutral 
IGM than faint LAEs, since the ionised spheres around these are 
typically larger. This indicates that the evolution of the LE happens 
mostly at the faint end, where a smaller neutral fraction is sufficient 
to prevent the observability. This scenario is similar to our observa¬ 
tions: most of the evolution between 2 = 5.7 and 2 = 6.6 happens 
at fainter luminosities, see Eig.[^ 

We use the following toy-model to show that a differential 
evolution indeed follows from some basic assumptions about how 
the process of re-ionisation happened. We first assume that the in¬ 
trinsic Lya LF at 2 = 6.6 (without IGM absorption) is similar to 
the 2 = 5.7 LF from |Ouchi et al.|j2008^ . Using this LF, we gener¬ 
ate a sample of LAEs with a minimum luminosity 10"^^ ® erg s“^ 
and consider three cases: First, we place them at random positions 
in our survey area, but also vary this to clustered positions, where 
the clustering has a functional form varying from Tow’ clustering 
(powerlaw) to ‘high’ clustering (exponential). In case of the pow- 
erlaw clustering, the probability p{x) that a LAE has a neighbour¬ 
ing LAE within a distance x, scales with p{x) ~ —x, while this 
scales with p{x) ~ e~^ in the case of exponential clustering. This 
means that for the powerlaw clustering, the expected number of 
neighbours within 1 Mpc is 2 times higher than based on a random 
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spatial distribution, while it is 13 times higher for the exponential 
clustering. 

Assuming that the LAEs are the only source of an ionised 
region around them (without any contribution from fainter neigh¬ 
bouring sources), we use the following formula for the radius of 
their Strdmgren sphere at 2 ; = 6.56 from|Cen & Haiman|(|2000|: 


D ,0 —1/3 f fesc.ion^j,em \l/3 / t 

(l.3xl05Ts-i) ( 2 ^ 77 ^' 


Mpc (6) 


where xigm is the neutral hydrogen fraction in the IGM, fesc,ion. 
the escape fraction of ionising radiation, N.y,em the number of ion¬ 
ising photons per second and t the time the emitter has been emit¬ 
ting ionising radiation. Under basic assumptions, we calculate the 
number of ionising photons by converting the Lya to Ha luminos¬ 
ity with: 


L{Ha) = 8.6fese,LyaL{Lya) erg s ^ (7) 

Here, f esc,Lya is the intrinsic Lya escape fraction and the 8.6 value 
corresponds to case B recombination ( |Osterbrock| 1989| >. Under the 
same case B recombination assumption and using a Salpeter IMF, 
the Ha luminosity is then converted to the number of ionising pho¬ 
tons using (e.g. lOrsi et al.|20T4) : 




L{Ha) 


1.37 X IO¬ 


CS) 


According to |Cen & Haiman| ( |2000^ , it takes an ionised 
sphere of roughly 1 Mpc for Lya to redshift out of resonance 
wavelength. Therefore, we use the above formalism to compute 
the ionised regions around LAEs, which can overlap if the sources 
are clustered, and check which LAEs can be observed, by being 
able to travel at least 1 Mpc through an ionised region through 
any sightline. As an output, we compile the observed Lya LF 
and compare this to the input. Our results depend obviously on 
our choice for the Lya escape fraction, the age of LAEs and the 
neutral fraction. However, we are only interested in qualitatively 
investigating the effect of this toy-model on the observed evolution 
of the LF, so the specific values are not that important, also because 
they all have the same powerlaw-scaling in Eq. 6. Our simulation 
likely breaks down when there is a strong dependence of intrinsic 
Lya luminosity on age, escape fraction or clustering (halo mass), 
since this would affect the relative sizes of ionised spheres and 
thus the observability of LAEs. 


After running our simulation 2000 times (to overcome 
sampling errors) for each set of clustering strengths and other 
input parameters, we find that the effect is that without clustering, 
brighter LAEs are more likely to be observed, since they are 
powerful enough to ionise their own surroundings. Changing any 
of the parameters except for clustering only varies the intrinsic 
Lya luminosity at which a LAE is able to ionise its surrounding 
enough. A stronger clustering leads to an increase in the number of 
observed faint LAEs, since they likely reside in ionised regions of 
larger LAEs, or are able to ionise a large enough region with some 
close neighbours themselves. 


Our results of the simulation are illustrated in the right panel 
of Fig. In this Figure, we fixed the parameters t, fesc,ion and 
fesc,Lyc to 100 Myr, 5 % and 30 %, respectively, and also fixed 
the clustering (to the maximum, exponential clustering strength) 
and only vary the neutral fraction. Our input LF (which is the same 


in the fully ionised case) is the LF at 2 = 5.7 from |Ouchi et al.| 
l |2008t . One can see that once the neutral fraction increases, the 
majority of the evolution happens at luminosities below 10"^® erg 
s“^, although the exact position depends on the specific neutral 
medium. This turnover point also depends on our choice of input 
parameters and therefore, since these are all uncertain, we do not 
use our model for estimates of the neutral fraction at 2 = 6.6. 
The turnover point corresponds to the luminosity at which a LAE 
is able to create a Strdmgren sphere of 1 Mpc. However, even in 
the most neutral fraction, there are still faint LAEs which can be 
observed. This is due to clustering, since they are located in the 
ionised regions of larger LAEs. The faint-end slope depends on 
the clustering, being the strongest if ionising sources are highly 
clustered. Concluding, the right panel of Fig. shows that in our 
simple toy-model for how re-ionisation happens and affects the 
observability of LAEs, the differential evolution that we observe 
between 2 = 5.7 and 2 = 6.6 can be explained by a higher neutral 
fraction of the IGM. It also shows that towards higher redshifts, 
where the Universe keeps becoming more neutral, strategies 
aiming at detecting LAEs benefit more from a wide-field approach, 
since LAEs are easier to be observed as they are able to ionise their 
own surroundings sufficiently. 

Note that our specific values for the neutral fraction are ar¬ 
bitrary since our input parameters are largely undetermined. We 
therefore do not aim to use our model to quantitatively measure 
the neutral fraction, but to qualitatively explain our observations. 
While typical models require much larger changes in the neutral 
fraction to significantly reduce Lya emission (e.g. |Dijkstra|2014t , 
our values are particularly low. This is because we assume that all 
Lya emission is absorbed when the ionised sphere is < 1 Mpc, 
while in reality only the blue part of the line-emission is absorbed, 
and the red wing broadens due to resonant scattering. We also ig¬ 
nore other processes which can increase the observability of Lya, 
such as galaxy outflows. 

The observed differential evolution for Lya selected galaxies 
is also seen in the spectroscopic follow-up of Lyman-break galax¬ 
ies. For example, |Ono et ak] ( |2012| ) finds that the drop in the frac¬ 
tion of UV selected galaxies which are detected in Lya at 2 > 6 
is stronger for UV faint galaxies, using a compilation of surveys. 
This is particularly interesting as at lower redshift, UV faint galax¬ 
ies tend to have a higher Lya fraction (e.g. |Stark et al.|2010^ . As 
mentioned by |Ono et al.H2012) , this differential drop indicates that 
re-ionisation happens first in over-dense regions and is qualitatively 
consistent with our observations of the Lya LF. 


6.3 Future surveys 

We can use our number counts, shown in Fig.|^ to estimate the ex¬ 
pected number of sources in future surveys such as the extragalac- 
tic Hyper Suprime-Cam (HSC) survey, which the Subaru telescope 
will undertake. Due to its large field of view (1.5 deg^), the HSC 
is suited excellently for pushing to a survey with a similar depth as 
the narrow-band data in UDS and a > 10 deg^ area, and even wider 
areas for slightly shallower surveys. From our spectroscopic con¬ 
firmed sample only (and without the filter correction), we find that 
the observed number density at 5 x 10"^^erg s“^, the luminosity of 
our brightest source in COSMOS, is roughly 30 times higher than 
based on the LF from |Ouchi et al.|p010| l (see Fig.[^. The HSC sur- 
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ve)[^has a planned deep component of ~ 30 deg^ (part of which is 
SA22) to a NB921 depth of 25.6 and an ultra-deep component (in 
UDS and COSMOS) of ~ 3.6 deg^ to a depth of 26.2. If we use our 
powerlaw fit from §5.1, we therefore expect that ~ 60 bright LAEs 
(> 5 X 10"^^erg s“^) will be found in total. These bright sources 
will be excellent targets for follow-up with JWST, to study the Lya 
escape fraction directly with Ha measurements, and to study the 
metallicity and ionisation state with other nebular emission lines. 
The HSC survey will also be able to study our brightest LAEs to 
more detail, because it will be able to observe the more extended, 
lower surface brightness regions because of deeper survey limits. 
Eurthermore, the HSC survey will also obtain additional wide ob¬ 
servations with the z = 5.7 Lya filter, such that the SA22 results 
at the highest luminosities can be compared. 

Another survey which recently started observations is the 
Javalambre Ph ysics of the Accele rating Universe Astronomical 
Survey (J-PAS) < Benitez et al. |2014| , which will survey 8500 deg^ 
with 54 narrow-band filters, from 3500 A up to a wavelength of 
10000 A, to a depth of ~ 22.5. Its goal is to study dark energy 
at redshifts z < 1, but it is also interesting for extragalactic stud¬ 
ies of quasars and LAEs due to its very wide area and large set 
of narrow-band filters. If we use our powerlaw fit (Table 7) and ex¬ 
trapolate it to a luminosity of 10^® erg s“^ (e.g. the Lya luminosity 
of the 2 = 7.085 quasar jMortlock et al.|2011) ), we estimate that 
it will find 230 bright LAEs per NB filter. A caveat however is that 
excluding interlopers is challenging with the depth of the optical 
imaging and without the availability of near-infrared. 

The brightest LAEs detected with optical surveys (e.g. up to 
2 ~ 7.3) are the best comparison with samples of even higher red- 
shift LAEs. To date however, no LAE has been detected in the near- 
infrared (e.g. 2 = 7.7 and 2 = 8.8; |Eaisst et al.| |2014[ |Matthee| 
|et al. 120141 1, but upcoming surveys will increase both the sensitivi¬ 
ties (which is a main difficulty due to the high sky background) and 
probed volumes. Nevertheless, our results, showing little to no evo¬ 
lution at the bright end, encourage further searches for the brightest 
LAEs, at least at 2 ~ 7 — 8. Eventually, the Euclid satellite will 
perform a deep spectroscopic survey over a wide area (40 deg^, 
5 X 10“^^erg cm“^ s”'^) in the near-infrared, which will detect 
bright 2 > 7.3 LAEs if they exist. 


are not as rare as previously thought, and have number densities of 
~ l.Slo.g X 10"® Mpc”®. 

• In our wide-field shallow SA22 data-set, we find 19 LAEs 
which are even brighter than the ones in UDS and COSMOS, al¬ 
though these have no spectroscopic follow-up yet. Even if there is 
significant contamination from low redshift interlopers or variable 
objects, all our results still hold. 

• After reproducing the analysis from |Ouchi et al.| j2010^ , we 
find that the main difference with our results is that we can more 
accurately include the brightest objects in our fit of the luminosity 
function and that we do a correction for our observational biases 
originating from the filter profile. We have additionally varied our 
EW criteria and used near-infrared data to exclude lower redshift 
interlopers, but find that this does not cause significant differences 
with the previous results for the ultra-deep surveys. Eor the deep 
COSMOS and shallower SA22 data, near-infrared data is however 
crucial. 

• Comparing our 2 = 6.6 data to lower redshift data at 2 = 
5.7, we confirm that there is evolution of the number density of the 
fainter LAEs. However, there is no evolution at the bright end and 
there is also evidence that the Lya LE deviates from a Schechter 
function at the brightest luminosities. 

• Based on our toy-model, we argue that this differential evolu¬ 
tion can be a sign of re-ionisation not being completed yet. Our 
model assumes that LAEs are only observed if they are in an 
ionised sphere which is large enough for the Lya photons to red¬ 
shift out of resonance wavelength, and that they are the only source 
of ionising radiation. Because of this, we preferentially observe the 
brightest LAEs, which have been able to ionise their own surround¬ 
ing enough to be observable. Eaint LAEs can only be observed if 
they are in the ionised spheres of more luminous LAEs, or when 
they are strongly clustered. 

• Additional wide field data at 2 = 5.7 and 2 = 7.3 will provide 
useful constraints on the evolution of the brightest LAEs around the 
re-ionisation epoch. 

• Finally, we use our results to make predictions for the upcom¬ 
ing Hyper Suprime Cam surveys at the Subaru telescope and also 
for upcoming very wide, cosmological surveys such as J-PAS and 
Euclid. 


1 CONCLUSIONS 

The main conclusions of this work are: 

• Using a combination of wide and ultra-deep surveys in SA22, 
UDS and COSMOS with the NB921 filter on the Subaru telescope, 
we obtain a large sample of LAEs at 2 = 6.6, spanning the largest 
dynamical range of luminosities to date, and derive a new luminos¬ 
ity function, which overcomes cosmic variance at the bright end 
because of our large and independent volumes. 

• We identify lower redshift interlopers (including extreme 
emission line galaxies at 2 = 1 — 2) with near-infrared data. This 
way we find a number of extreme emission line galaxies with rest- 
frame EWs of > 400 A. 

• In total, we find 99-l-15-tl9 (UDS, COSMOS, SA22) LAE can¬ 
didates at 2 = 6.6. Of these, 18 are confirmed spectroscopically. 16 
in UDS bylOuchi et al.|(|Mol and 2 in COSMOS by|Sobral et al.| 
(|20T^. These two sources are the brightest LAEs at 2 = 6.6 so 
far and confirm that sources with a luminosity similar to Himiko 

4 www.naoj.org/Projects/HSC/surveyplan.html 
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